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SUMMARY 

Numerical solutions of the unsteady Euler equations are obtained using the 
classical fourth order Runge Kutta time marching scheme. This method Is fully 
explicit and Is applied to the governing equations In the finite volume, con- 
servation law form. In order to determine the efficiency of this scheme for 
solving turbomachinery flows, steady blade-to-blade solutions are obtained for 
compressor and turbine cascades under subsonic and transonic flow conditions. 
Computed results are compared with other numerical methods and cascade tunnel 
measurements. The present study also focuses on other Important numerical 
aspects Influencing the performance of the algorithm and the solution accuracy 
such as grid types, boundary conditions, and artificial viscosity. For this 
purpose, H, 0, and C type computational grids as well as characteristic and 
extrapolation type boundary conditions are Included In the solution procedure. 


INTRODUCTION 

A broad spectrum of flow conditions are encountered at various stages of 
a modern gas turbine engine operational cycle. These conditions are highly 
complex and are significantly diversified from one component to the other. 
Development of analytical tools closely representing the Important flow fea- 
tures Is a necessary step for efficient design, operation, and performance 
Improvements. Due to the lack of computing resources and limitations of the 
existing analytical methods, no single study can hope to provide all the Infor- 
mation and answers to the complete needs and satisfaction of the design engi- 
neer. Hence, development of numerical methods Is customarily carried out In 
stages, ultimately hoping to produce a model representing the gas turbine flow 
field as realistically as possible. 

The purpose of the present study Is to develop a fast and accurate Euler 
cascade computer flow code that can be used for the design and performance 
Improvements of critical gas turbine engine components such as compressors and 
turbines. Development of the solution procedure Is based on the Jameson's 
fourth order Runge Kutta numerical Integration scheme (ref. 1). This scheme 
has been successfully applied and shown to yield fast results with good accu- 
racy for external aerodynamic flow problems. Other major advantages that 
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prompted the authors to apply this method for Internal turbomachinery flows are 
the simplicity and straightforward Implementation procedure. Independence of 
the steady-state solution on the time step taken and the vectorlzable nature 
of the scheme. 

The numerical results obtained In this study address several Important 
Issues related to the performance of the present method compared to the other 
schemes that exist In the literature for cascade flows (refs. 2 to 8). These 
Issues Include solution convergence, mass and total pressure conservation, grid 
types such as H, 0, and C, boundary condition methods, and Influence of flow 
geometry variation such as compressors and turbines. 


GOVERNING EQUATIONS AND NUMERICAL INTEGRATION SCHEME 


The Euler equations governing the two-dimensional Invlscld flows can be 
written In the Integral form as 
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where AR Is the flow domain, R Is the flow boundary, x and y are the 
Cartesian coordinates. The vectors U, F and G In equation (1) are 
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Here, p, u, v, p, E, and H are the density, velocity components In the x 
and y directions, pressure, total energy, and total enthalpy respectively. 
For a perfect gas 
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and 


where y Is the specific heat ratio. Approximation to equation (1) leads to 
finite volume numerical method In conservation form. 


Runge Kutta Integration Scheme 

The flow domain Is divided Into a number of small finite volume (finite 
area In two dimensions) computational cells and the numerical approximation of 
equation (1) Is applied to each cell separately. This procedure leads to a 
system of ordinary differential equations of the form 
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^ [U • AA] + Q*U - DU = 0 


(5) 


which can be solved by any number of Integration schemes. In this equation, AA 
Is the cell area, Q Is the spatial differencing operator and D(U) Is the 
added artificial dissipative term to supress numerical Instabilities In the 
region of very large flow variations. The exact form of equation (5) when 
applied to the x momentum component for example. Is (with the dissipative 
term omitted for simplicity) 
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where the flux velocity Qj< Is 
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and the sum Is over all four sides of the cell. In figure 1 the four sides of 
the cell (1,J) are denoted by numbers one to four. The dependent flow varia- 
bles such as p, pu, pv, etc. are cell centered and the values ( pU) ^ across 
any side k are average values on the two sides of the face. 

For example, 

<pU >l = j[(P"),,J - ( p “)l,j.i] 

(pu) 4 > j [ ( ' >u) u * (pU) 1.1 J ] 

AX 1 = (x D - x A ) . T ; Ay 4 = (y c - y Q ) • J (8) 


and so on. 


Equation (5) Is Integrated using the modified four stage Runge Kutta 
scheme In which the dissipative terms are frozen at the values of the first 
stage. The values at the time level "n" are then updated to the new time 
level "n ♦ 1" In the following four stages: 
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where At Is the local time marching step. The dissipative operator D Is a 
blend of second and fourth differences. The operator depends on the local 
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flow behavior and Is Included In the exact fashion described by Jameson 
(ref. 1). 


BOUNDARY CONDITIONS 

With reference to figure 2, there are four types of boundaries that con- 
stitute the flow region and require special treatment. These are the Inflow 
boundary AH, the outflow boundary DE, the solid-wall boundaries (BC and 
GF), and the periodic boundaries AB, CD, HG, and FE. The appropriate set of 
conditions for each of these boundaries Is described below. 

Inflow and outflow boundary conditions used In this study were the char- 
acteristic and extra-polatlon type. If the flow Is subsonic there will be 
three Incoming (right running) characteristics and one outgoing (left running) 
characteristic at the Inflow while the opposite Is true at the outflow boundary 
where there are three outgoing (right running) characteristics and one Incoming 
(left running) characteristic. By the theory of characteristics, three condi- 
tions may therefore be specified at the Inflow and one condition at the outflow 
(ref. 9). The remaining conditions are numerically determined by the solution 
of the differential equations. The three conditions specified at the Inflow 
are the total pressure, total temperature, and the flow angle. The remaining 
condition Is obtained by extrapolating the outgoing or upstream running Rlemann 
Invariant from the Interior to the Inlet. For extrapolation type, the static 
pressure Is extrapolated from the Interior cell to the Inflow. For supersonic 
axial Inflow, the Inlet Mach number Is also specified and held constant In 
addition to the three conditions mentioned above. 

At the outflow boundary, the one physical condition specified Is the 
static pressure corresponding to the desired exit Mach number. The three 
numerical conditions come from extrapolating the downstream running Rlemann 
Invarlent, the y-veloclty component v and the total energy E. For super- 
sonic axial outflow, the static exit pressure Is also extrapolated from the 
Interior point. 

On the blade surface, the "zero flux" conditions are Imposed. However, 
for estimating the contributions to the momentum equations due to the pressure 
terms, we need to evaluate the pressure on the blade surfaces. The normal 
pressure gradient at the boundary cell centers Is determined by the same fash- 
ion described In reference 1. Extrapolation using the pressure gradient and 
the pressure at the cell center determines the pressure at the wall. 

The periodicity Is Imposed by setting the flow variables In the cell cen- 
ters that lie within the computational domain equal to the corresponding per- 
iodic cells that lie outside the computational domain. For example (with 
reference to fig. 2), the unknown flow quantities In the control volumes that 
lie Immediately above the periodic boundaries AB and CD are set equal to the 
corresponding cell values that lte Immediately above the periodic boundaries 
HG and FE whose values are known from the Interior point calculations. 
Additional details on the Implementation of these boundary conditions can be 
found In references. 1 and 7. 
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RESULTS AND DISCUSSION 


Numerical solutions of the Euler equations using the Runge Kutta scheme 
described earlier were obtained for selected cascade test cases with widely 
varying geometries and flow conditions. Results are presented and compared 
with solutions of other numerical methods and experimental data. Computational 
grids of different types were also used In the calculations to determine their 
Influence on the solution accuracy. 

The Mach number contours calculated for a transonic flow through a com- 
pressor cascade are Illustrated In figure 3. The number of grid points used 
In the computations are 65 by 15. The cascade geometry and the grids are shown 
In figure 2. The flow Is subsonic at the Inlet and accelerates to a peak Mach 
number of 1.246 on the blade suction side before shocking down to an exit Mach 
number of 0.67. The computed surface Mach numbers agree very well with the 
predictions of other numerical methods. However, there are no experimental 
data available for this case. 

Figure 4 depicts a high work turbine guide vane geometry (ref. 13) and the 
80 by 17 computational grids used for obtaining a fully subsonic flow through 
the passage. The computed Mach number contours are shown In figure 5. It Is 
seen that the flow accelerates from an Inflow Mach number of 0.11 to an exit 
value of 0.84 In a very short distance. The surface static pressure ratios 
obtained by the present method are shown In figure 6 along with available test 
data (ref. 13) an another code (ref. 14) using the two step explicit MacCormack 
numerical scheme. Considering the type and coarseness of the grids employed 
for this difficult test case, the present code yielded results with very good 
accuracy. Better accuracy can be achieved by refining and Increasing the num- 
ber of grid points In the flow field. 

Next example Is a flow through a NASA turbine stator rows shown In 
figure 7. The calculated critical surface Mach numbers are shown In figure 8. 
The wind tunnel data (ref. 10) and the numerical results of another code 
(ref. 3) computed In house are also shown In this figure for comparison. The 
agreement of the present results with the test data as well as with the other 
method Is very good. However, the poor agreement of the two numerical results 
with the test data near the leading edge Is mainly due to the nature of the 
grid In this region. The numerical schemes treat the leading edge with "cusps" 
for solution convergence which Is not present In the actual airfoil. 

The treatment of leading and trailing edges with artificial cusps for 
numerical convenience Is a common feature while computatlng flows with sheared 
"H" type grids as done In the previous test cases. This Is particularly true 
for cascades with high Incidence, large flow turning, and thick leading/ 
trailing edge profiles. Numerical results obtained with "H" grids are often 
Inaccurate and misleading In these Important flow regions. For Improved accu- 
racy, It Is preferred that computations be performed on body fitted orthogonal 
or near orthogonal "C" or "0" type grids (ref. 11). To Illustrate this 
point, the surface Mach numbers calculated on a "C" grid for the same test 
case of figure 7 are compared with the test data and this Is shown In figure 9. 
It can be seen that the agreement Is excellent and far Improved near the lead- 
ing and trailing edges compared to the "H" grid results. 

Figure 10 shows an AVC0 rotor blade and the "C" type computational mesh 
for which numerical solutions were obtained. This case was particularly chosen 
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to test the ability of the present code to calculate transonic flows through 
turbine blade passages with high flow Incidence, large turning, and thick blade 
profiles at the leading and trailing edges. The Inlet flow angle Is 55°, the 
Inflow Mach number Is 0.66 and the outflow Mach number Is near transonic. The 
Mach number contours near the leading edge region Is shown In figure 11 to 
Illustrate the quality of results that can be obtained on a "C" grid with the 
present method. Computed surface Mach numbers are plotted In figure 12 along 
with the numerical results of reference 3. The agreement between the two 
methods are very good except near the shock location where the present method 
predicts a peak value of 1.108 compared with 1.152 obtained by the Denton's 
code (ref. 3). 

The VKI gas turbine rotor blade and the computational grids on which tran- 
sonic flow solutions were obtained Is pictured In figure 13. Figure 14 com- 
pares the predictions of surface Mach numbers with test data obtained from two 
different cascade tunnels and two different exit flow conditions. The agree- 
ment between the present calculations and the two sets of test data Is excel- 
lent especially for the case corresponding to exit Mach number of 1.19. 


CONCLUSIONS 

A computer code for solving the Euler equations using the four stage Runge 
Kutta Integration scheme has been developed for turbomachinery flow field cal- 
culation. The program has been successfully applied to predict blade-to-blade 
flows for many cascades with different geometries and flow conditions. Numer- 
ical results Indicate that the present method can be applied to yield fast 
results with good accuracy for a wide variety of cascade configurations and 
flow conditions. The "C" type grids produce the best overall results for any 
particular test case from the stand point of accuracy, simplicity of Implemen- 
tation, and boundary condition treatments. The present code, which could be 
used In conjunction with any type of user opted computational grids, Is simple, 
efficient, and accurate enough to be used for cost effective preliminary and 
advanced aerodynamic design studies. 
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Figure 2. - Compressor cascade and computational 
grids. 




Figure 4. - Turbine guide vane and computational mesh. 
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Figure 6. - Calculated surface static pressure vari- 
ation. 
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Figure 12. - Surface Mach number variation for AVCO 
cascade. 



Figure 13. - VKI turbine cascade and computational 'C' grids. 
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Figure 14. - Surface Mach number distribution. 
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